clear all;
close all;




a1 = load('LinRegpred.txt');
HR=a1(:,2);
HRE=a1(:,3);
X=a1(1000:1500,2:3);
HRN=[HR(1:1000) HRE(1:1000)];

step=30

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% phase 1: lower with lower estimation
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for i=1:step
    HRA(i,1)= randi([0,60],1,1);
   
end

for i=1:step
    HRA1(i,1)=  HRA(i,1) + sign(rand(1,1))*randi([0 30],1,1); 
end

HRAS1=[HRA HRA1];

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% phase 2: higher values with higher estimation
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear HRA HRA1

for i=1:step
    HRA(i,1)= randi([125,220],1,1);
   
end

for  i=1:step
    HRA1(i,1)=  HRA(i,1) + sign(rand(1,1))*randi([0 30],1,1); 
end

for  i=step:2*step
    HRA(i,1)= randi([125,220],1,1);
   
end

for   i=step:2*step
    HRA1(i,1)=  HRA(i,1) - randi([20 50],1,1); 
end

HRAS1=[HRAS1; HRA HRA1];

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% phase 3: normal values with lower estimation
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear HRA HRA1

for  i=1:step
    HRA(i,1)= randi([60,120],1,1);
   
end

for i=1:step
    HRA1(i,1)=  HRA(i,1) - randi([20 60],1,1); 
end

HRAS1=[HRAS1; HRA HRA1];

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% phase 4: normal values with higher estimation
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear HRA HRA1

for  i=1:step
    HRA(i,1)= randi([60,120],1,1);
   
end

for  i=1:step
    HRA1(i,1)=  HRA(i,1) + randi([20 60],1,1); 
end

HRAS1=[HRAS1; HRA HRA1];

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% phase 5: normal values with normal estimation
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


clear HRA HRA1

for  i=1:step
    HRA(i,1)= randi([60,120],1,1);
   
end


for  i=1:step
    HRA1(i,1)=  HRA(i,1) + sign(rand(1,1))* randi([0 10],1,1); 
end


m=1;
for  i=1:step
  moi(m,1)=  HRA(i,1);
  moi(m,2)= HRA1(i,1);
  m=m+1;
end

HRAS1=[HRAS1; HRA HRA1];
clear HRA
HRA=HRAS1;

[l,c]=size(HRA);
A=-1;
B=1;
HRA(:,3)=A;
HRA(151:end,3)=B;




figure(20), plot(HRA(:,1),HRA(:,2),'ro')
 hold on
 plot(moi(:,1),moi(:,2),'b*')
 
%  figure(21)
%  plot(HRA(1,A==-1),HRA(2,A==-1),'ro',...
%     HRA(1,A==1),HRA(2,A==1),'b*')
return







for i=1:500
    HRA(i,1)= randi([10,60],1,1);
   
end

for i=500:1000
    HRA(i,1)= randi([120,200],1,1);
end


for i=1:1000
    HRA1(i,1)= HRA(i,1) + sign(rand(1,1))*randi([20 100],1,1); 
end
%A=[HRM; HRE];

HRA=[HRA HRA1];

figure(20), plot(HRA(:,1),HRA(:,2),'ro')
return


plot(X10(bel==1,1),X10(bel==1,2),'ro',...
    X10(bel==2,1),X10(bel==2,2),'b*')



return 
mv_ini=[mean(HRN); mean(HRA)];
mv_ini=mv_ini';
X=[HRN; HRA]
cA=cov(HRN)
cN=cov(HRA)
mc_ini(:,:,1)=cA';
mc_ini(:,:,2)=cN';
e=0.001
maxiter=300
sed=0
[ap,cp,mv,mc, iter,diffvec] = GMDAS(X',mv_ini,mc_ini, e,maxiter, sed)
[qw,bel]=max(cp')
figure(20),plot(X(bel==1,1),X(bel==1,2),'ro',...
    X(bel==2,1),X(bel==2,2),'b*')

X=a1(1000:1500,2:3);
X(490,1)= X(490,1) + 100;
for i=1:501
   d1= compute_gauss_dens_val(mv(:,1),mc(:,:,1),X(i,:)')
%    d2= compute_gauss_dens_val(mv(:,2),mc(:,:,2),X(i,:)')
%    d=[d1 d2]; 
%    [d,ind]=max(d)
%     dodo(i)=ind;
end


return

% Plot the data set
X10=HR(1:100,:);

clear all
N(:,6)=[];
A(:,6)=[];
A=A'
N=N'
mN=mean(N)
A=A'
N=N'
mN=mean(N)
mA=mean(A)
cA=cov(A)
cN=cov(N)
mv_ini=[mA' mN']
mc_ini(:,:,1)=cA'
mc_ini(:,:,2)=cN'
e=0.1
maxiter=300
sed=0
[ap,cp,mv,mc, iter,diffvec] = GMDAS(X',mv_ini,mc_ini, e,maxiter, sed)
X=[A; N]
[ap,cp,mv,mc, iter,diffvec] = GMDAS(X',mv_ini,mc_ini, e,maxiter, sed)
[qw,bel]=max(cp')
e=0.001
[ap,cp,mv,mc, iter,diffvec] = GMDAS(X',mv_ini,mc_ini, e,maxiter, sed)
[qw,bel]=max(cp')
Untitled
X=[ N;A]
[ap,cp,mv,mc, iter,diffvec] = GMDAS(X',mv_ini,mc_ini, e,maxiter, sed)
[qw,bel]=max(cp')
Untitled
A=A(:,3:4)
N=N(:,3:4)
figure(20),plot(X10(bel==1,1),X10(bel==1,2),'ro',...
    X10(bel==2,1),X10(bel==2,2),'b*')



return 



%return




% X10=X10'
% %X10=[X10'; X12']'
% %%plot(X10(1,:),X10(2,:),'k.')
 %X10(800:5:1000,:)=X10(800:5:1000,:)+ 100;

[bel,C]=kmeans(X10,2)
figure(20),plot(X10(bel==1,1),X10(bel==1,2),'ro',...
    X10(bel==2,1),X10(bel==2,2),'b*')

[A,C]=MinVolEllipse(X10', 0.01);
hold on
Ellipse_plot(A,C)
%ezplot((A-repmat(C,1,2))'*X10'*(A-repmat(C,1,2))=1)
%return 

%X11(990:5:1000,:)=X11(990:5:1000,:)+ 10;
%X12=diff(X11(:,1))- repmat(mean(X10'),size(X10,2),1)';
%X13=diff(X11(:,2));
%X10=[X12'; X13'];
%X10(2,100)=150;
%X10(2,200)=50;
%  medano=median(X10)
% 
% 
% v=std(X10);
% X10=X10- repmat(medano,1000,1);
% X10= X10 ./ repmat(v,1000,1);

[l,N]=size(X10);
figure(1), plot(X10(:,1),X10(:,2),'k.')
figure(1), axis equal
%return
clear bel
% 2. Apply the spectral clustering algorithm
e=2; %Thershold for the distance in the definition of W. Also try e=3
sigma2=2;%sqrt(10); %The sigma^2 in the exponential in the definition of W
bel=spectral_Ncut2(X10',e,sigma2);

% Plot the clustering results (see Figure 7.10(a))
figure(2),plot(X10(bel==0,1),X10(bel==0,2),'ro',...
    X10(bel==1,1),X10(bel==1,2),'b*')
figure(2), axis equal

 mean(X10(bel==0,:))
  mean(X10(bel==1,:))